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Thermoacoustic instability occurs in many combustion systems, such as aero-engine afterburners, rocket 
motors, ramjets and gas turbines. It most often arises due to the coupling between unsteady heat release 
and acoustic waves. In this work, nonorthogonality analysis of a choked combustor with a gutter confined 
is conducted. Such configuration is used as a simplified model of the afterburner of an aero-engine. A 
thermoacoustic model is developed first to study the nonnormal interaction between acoustic distur¬ 
bances and a premixed V-shaped flame anchored to the tip of the gutter. Eigenmode nonorthogonality 
analysis is then conducted. The thermoacoustic system is shown to be nonnormal and characterized 
by nonorthogonal eigenmodes. The nonorthogonality is identified to arise from both the complex bound¬ 
ary condition and the monopole-like flame. However, the contribution from the Robin-type boundary is 
approximately 1.5% of that from the flame. Thus the flame is identified to play a dominant role. One prac¬ 
tical conclusions is that acoustic disturbances undergo transient growth in a combustion system with 
nonorthogonal eigenmodes. Such finite-time growth, which cannot be predicted by using classical linear 
theory might trigger high-amplitude self-sustained oscillations. 

© 2014 Elsevier Ltd. All rights reserved. 


1. Introduction 

Thermoacoustic instabilities [1-3] are characterized by large- 
amplitude self-excited pressure oscillations [4-13] caused by a 
coupling between acoustic waves and unsteady combustion. The 
acoustic waves perturb the flame and change the instantaneous 
heat release rate [14,15], Since the unsteady heat release is an 
efficient monopole-like sound source [16], the acoustic waves 
can gain energy from their interaction with unsteady combustion. 
This feedback can result in the pressure oscillation amplitudes suc¬ 
cessively increasing i 17,18]. Eventually, some nonlinearity in the 
combustion system will limit the amplitude of oscillations [19], 
The combustion-excited oscillations might be so intense as to 
cause structural damage and costly mission failure 1], 

The occurrence of thermoacoustic instabilities has been a 
plaguing problem in the development of combustors for rocket 
motors [1], ramjet [20], the afterburner of aero-engines, and 
land-based gas turbines [1], Over the last decades, thermoacoustic 
instabilities have been the subject of intense research activity, aim¬ 
ing to better understand and predict them [21], Linear stability 
analysis [22] is generally conducted for studying the stability 
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behaviors of combustion systems via calculating eigenfrequencies. 
A combustion system is linearly stable, when all eigenmodes decay 
exponentially. However, when the eigenmodes are not orthogonal, 
the combustion system is nonnormal and acoustic disturbances 
present will undergo transient growth [23-26]. If the transient 
growth rate is large enough, thermoacoustic instability might be 
triggered by causing small-amplitude acoustic disturbances to 
grow to amplitudes high enough to make nonlinear effects 
significant. 

Nonnormality analysis of a solid rocket motor was attempted 
by Mariappan and Sujith [24]. However, the boundary conditions 
are modeled as acoustically closed at both ends for simplicity. This 
means that the rocket motor boundary conditions are Dirichlet- 
type ones. It was shown that transient growth of acoustic distur¬ 
bances occurred in such nonnormal system. Similarly, a Rijke-type 
combustion system with Dirichlet-type boundary conditions was 
investigated [23,25], A heat source is modeled and simply 
described by using time-lag formulation. The interaction between 
a flame and the acoustic waves was not discussed. This partially 
motivated the present work. 

In this work, nonorthogonality analysis of a longitudinal 
combustor with a choked inlet and an open exit are conducted. A 
V-shaped premixed flame anchored on a bluff-body/flame-holder 
is considered and it is assumed acoustically compact. The flow 
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Nomenclature 


Af flame surface area, m 2 

A coefficient matrix 

B coefficient matrix 

c u . c d speed of sound upstream and downstream of the flame, 
m/s 

c p , c„ heat capacity ratio at constant pressure and volume, 
kj/kg K 

Cp identity vector 

G G-equation 

h specific enthalpy, W/m 2 K 

H Heaviside step function 

j j = V 3 ! 

M u , M d Mach number upstream and downstream of the flame 
M interaction index 

n the unit vector of the flame propagation 

p' pressure fluctuation, Pa 

p mean pressure, Pa 

q heat flux, W/m 2 

Q instantaneous heat production per unit volume J/m 3 

Q', Q fluctuating and mean heat release rate, kj/s 

r radius of the flame front, m 

R gas constant, J/mol K 

R a , R b radius of the burner and the tube, m 

S t flame Strouhal number 

T Temperature, K 

T r flame transfer function 

Uf unburned gas velocity, m/s 

u u , u d mean velocity upstream and downstream of the flame, 

m/s 

Uf , (if mean and fluctuating part of velocity at flame position, 

m/s 

Vf flame speed, m/s 

V N state vector in Eq. (41 ) 

W coefficient matrix in Eq. (41) 

z axial position, m 

Zf heat source location along z, m 


z u , z d axial position of inlet and outlet, m 

Z acoustic impedance, m/s 

Z u , Z d reflection coefficient at inlet and outlet 

Greek letters 

X coefficient in Eq. (22) 

d(z) dirac delta function 

e coefficient in Eq. (15) 

y the ratio of specific heat 

kt weighting factor corresponding to the trapezoidal inte¬ 

gration 

</>o initial phase, rad 

m oscillation frequency, rad/s 

a m grow rate of mth eigenmode 

X Rayleigh index 

p air density, kg/m 3 

£(r, t) flame front location 

x v time delay between heat and pressure fluctuation, s 

r u , x d time delays, s 
Q m eigenvalues 

Re real part 

Im imaginary part 

Subscript 

u upstream 

d, downstream 

o branch outlet 

/ flame 

m mth eigenmode 

Superscript 
- mean value 

A Fourier transform 

* complex conjugate 

/ fluctuation part 


fluctuations are expanded using decomposed traveling plane 
waves. In Section 2, the thermoacoustic model equations are devel¬ 
oped. The flame front is tracked by using the classical G-equation. 
In Section 3 an analytical method to quantify eigenmodes nonor¬ 
thogonality is introduced. In Section 4, the flame response to 
oncoming acoustic disturbances are discussed in terms of the flame 
transfer function and the time evolution of the flame front. In addi¬ 
tion, the contributions to the eigenmodes nonorthogonality from 
the flame and the choked boundary are compared and discussed. 
To gain insight on the effect of Robin-type boundary condition 
on the nonorthogonality, the combustor in the absence of the flame 
is considered. 


2. Thermoacoustic model of the premixed choked combustor 

The numerical model considers a confined premixed V-shaped 
flame burning in the recirculation zone of a bluff body (flame- 
holder) in a longitudinal combustor. The combustor is associated 
with a choked inlet and an open exit. Such configuration is one 
of the simplest generic geometries of an aeroengine afterburner 
with Robin-type boundary condition [26,27], It has been exten¬ 
sively investigated experimentally [28,29]. However, nonorthogo¬ 
nality and nonnormality analyses of such combustion system 
have not been reported yet. For this, a thermoacoustic model with 
a V-shaped turbulent flame is developed. The numerical model 


captures two physical processes; the generation and propagation 
of the acoustic waves within the combustor, and the dynamic 
response of the flame to the acoustic waves. Here, the flame front 
is tracked in real-time by the kinematic G-equation, of which a 
time-invariant turbulent flame speed is used. 

2.1. Model for combustion-driven acoustic waves 

A premixed flame with constant fuel-air ratio burning in a 
cylindrical combustor is considered, as shown in Fig. 1. The geom¬ 
etry of the combustor and flow conditions are summarized in 
Table A.l in Appendix A. The flame, stabilized in the wake of an 
axisymmetric flame-holder (center-body), is acoustically compact 



Fig. 1. Schematic of a cylindrical duct with a premixed V-shaped flame anchored on 
the flame-holder/gutter. 
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and it is modeled as a thin sheet. The presence of the flame results 
in the mean temperature T, the speed of sound c and the gas den¬ 
sity p undergoing a sudden jump. Thus, we use subscript u and d to 
denote the parameters in the pre-burn region and after-flame 
region respectively. The length of the combustor is L and the flame 
is being axially placed at Z/. The pressure of the gas mixtures in the 
combustor is always low enough to justify the assumption of the 
perfect gas laws. Thus if p denotes the pressure, T the temperature 
and p the density, p = RTp = pTR u /M w . Here, M w is the mean 
molecular mass of the mixture and R u is the universal gas constant. 

Consider a control volume enclosing the acoustically compact 
flame, the density p varies through changes in both the pressure 
p and specific enthalpy h. Thus 

dp_ldp dp d h 
dt c 2 df dh p dt 

For a perfect gas p = p/RT = yp/c 2 , dp/dh\ p = -p(y - 1 )/c 2 and so 
Eq. (1) can be rewritten as 


dp _ 1 dp p(y - 1) d/i 
dt c 2 dt c 2 dt 


where y is the ratio of specific heat. When viscous and heat conduc¬ 
tion effects are neglected, pd/i/dt = Q, where Q(t) = Q + Q'(t) is the 
instantaneous heat input/unit volume. An over-bar denotes a mean 
value and a prime denotes a perturbation. 

The mean flow has vanishingly small Mach number (i.e. 
M 2 «l .0), so that the flow properties such as density p and pres¬ 
sure p can be assumed to consist of a mean and a fluctuating part 
as p = p + p' and p = p + p'. The mass and momentum equations 
over the control volume, when linearized with respect to the per¬ 
turbation quantities, can be written as 


1 dp' _ du' 
c 2 dt + ^ dz 




(3) 


_ du 1 dp' 

fut+is - 0 


(4) 


By differentiating Eq. (3) with respect to time t, Eq. (4) with respect 
to z and subtracting, a general combustion-driven wave equation is 
obtained namely, 


1 c>V d 2 p' 

c 2 dt 2 dz 2 


(z~zf)) 


(5) 


where S(z) is dirac delta function describing localized compact heat 
source. When there is no heat input, i.e. Q(t) = 0, Eq. (5) reduces to 
the well-known acoustic wave equation: 


1 d 2 p' d 2 p' 
c 2 dt 2 dz 2 


( 6 ) 


If the acoustic waves in different combustion regions are modeled 
as right- and left traveling waves p + (z, t) and p~(z, t), as shown in 
Fig. 1, then the instantaneous pressure fluctuation upstream and 
downstream of the flame can be described as 


where M u and M d denote the Mach numbers of upstream and 
downstream of the flame. Manipulating Eqs. (7) and (8) gives rise 
to the acoustic traveling waves as given as p ± (z,t) = (p'(z,t)± 
pcu'(z, t))/2. 

The combustor boundaries are modeled using pressure reflec¬ 
tion coefficients [25], At the open downstream end, the reflection 
coefficient Z d (co) is used to couple the reflected wave p d (z d , t) with 
the incident one pf(z d ,t) as given in Eq. (9a). And at the choked 
upstream end, Z u (co) is used to couple p+(z d , t) to p u (z d . t ), as given 
in Eq. (9b). 

p d (z d , co) = Z d {co)p+ exp (- 2jm(Zd — (9a) 

V ^0 ~ M d)J 


p+(-z u ,co) = Z u (co)p u exp 


2}co(Zf - z u ) \ 

MI-MI)) 


(9b) 


Since the mass flow rate at the choked inlet is constant, it can be 
then shown that the reflection coefficient Z u = 1 - M„/l + M u . 
The downstream reflection coefficient Z d is set to be to -1.0 for 
simplicity. 


2.2. Model for premixed V-shaped flame dynamics 


The inhomogeneous combustion-driven acoustic wave equation 
as given in Eq. (5) can be solved if Q'(t) is known. Following the 
work [27], the unsteady heat release Q'(t) from the premixed 
ducted flame imitates the evolution of the flame surface area ratio 
such that 


Q.\t) __ A}(t-T) 

Q * A f 


( 10 ) 


To capture the dynamics of the V-shaped flame response with 
oncoming acoustic disturbances, the classical G-equation is used 
for flame-front tracking [ 27 ] to study the nonnormal interaction 
between the flame and flow disturbances. The main features of 
the flame model are reproduced here for completeness and integ¬ 
rity. The flow is assumed to be axisymmetric and combustion is 
assumed to occur on a thin surface whose axial position at radius 
r is given by z s = {( r.t ) as shown in Fig. 1. The flame front corre¬ 
sponding to a particular level G(Z/, r.t) = 0 of a scalar field G 
describes the flame initiation surface, as, Z/ = {(r, t). The flame sur¬ 
face is assumed to propagate itself in the direction of its normal 
n/ = VG/|VG|, into the unburnt fluid at a constant burning velocity 
Vj relative to the unburnt gas. That is to say the surface 
G(Zf,r,t) = 0 moves in the direction of its normal n/ with speed 
u f Hf - V f , where U/ = (U/, vf) is the unburned gas velocity, and 
n/ points downstream as shown in Fig. 1. Mathematically this state¬ 
ment is equivalent to AG/ df = 0: 

Tt = Tt + {u f~ V f n/) ' VG (11) 


®fsii=! p+ [ p “( t_ «) +Pu ( t+ «)] Zu ^ 2/ 

z r <z<z * 

The instantaneous velocity fluctuation is given as 

u{z . t) = {“ u+ ^[ Pj ( t_ ^)" P “( t + ^fc)] Zu<Z<Z/ 

z ' <z<Zd 


(7) 


( 8 ) 


After substituting n f = VG/|VG| and G(Zf.r.t)=Zf-g(r,t), the 
familiar G-equation is obtained, which describes the location of 
the flame surface [27]: 


di(r,t) 

dt 


+ v f 


dj(r,t) 

dr 


+ V f 


1 + 


m 



= Uf 


( 12 ) 


With the assumption that the density change across the surface 
G(Z/, r, t) = 0 is negligible [27], the particle velocity is given by that 
in the oncoming flow (u/,0) and the equation describing the time 
evolution of {(r, t ) is shown as 
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dj(r,t) 

dt 


= u t (Zf, t) - Vf 


1 + 


V dr 


1/2 


(13) 


It is apparent that the flame front evolution is dynamically coupled 
with the unsteady velocity u/(t) at the flame-holder, which is 
related to the oncoming acoustic disturbance velocity u u as 


u f (t) = 


K 

(Rl-R 2 a) L 


-^r(Pt(Zf,t)-Pu&,t)) 


(14) 


where R a and R b are the radius of the flame-holder and the combus¬ 
tor respectively. When u/(t) > Vf, the flame is attached at the 
flame-holder, i.e. £(R a , t) = 0. And Eq. (13) can be rearranged to give 


dj(r,t) 

dr 


uj(t) 


\ Vj 


- 1 on r = R„ 


(15) 


However, when u f (t) < V f , the flame no longer remains to attach at 
the flame-holder and it propagates upstream. Thus the boundary 
condition on r = R a is given as d((r,t)/dr= 0, if u s (t)<V s or 

t(Ra t) < 0. 

In our calculation, £(r, t) is evaluated at 300 locations across the 
combustor radius and a 4th order Runge-Kutta algorithm is used 
for the numerical time integration. 

Once £(r, t) is known, the instantaneous flame surface area A f (t) 
can be readily calculated by 


Af(t) : 


r 2n rR b 

dd r| 1 + 
Jr„ 


d £(r, t) 
dr 


1/2 


dr 


(16) 


Since ((r,t) = £(r) + {'(r, t), substituting it into Eq. (16) and using 
binomial expansion scheme leads to 


w ‘i>/,T + ( 

-s: 


. 2 — — + (—^ dr 
dr dr + \dr ) 


2nr\j 1 




I)’- 


r R t 

L 


2nr 


K_ d£ 
dr dr 


--dr 


(17) 


It can be seen from Eq. (17) that the mean and linearized fluctuating 
values of the flame surface area are given as 


Af = 


j 

JR a 

f R i 

w-L 


2nr\ 1 


2nr 




e?(r,t) 9|(r) 

tlr . a L . dr 


1 


{¥) 


(18) 


(19) 


The mean value of Eq. (13) leads directly to u f = V f <J 1 + dc(r)/dr. 
Manipulating it leads to dc(r)/dr= sjiij - V 2 f /Vf. After integration 
and application of boundary condition, it can be shown that 


- 1 


i(r) = (r-R a ) ] 

Substituting Eq. (20) into Eqs. (18) and (19) leads to 

_ n(R 2 b -R 2 a )u f 

A f = —-- 


( 20 ) 


( 21 ) 


where / = -2n^J\ - Vj/uj. To simplified Eq. (22), integration by 
parts and the boundary conditions are applied. If the flame front 
is discretized into K flame elements equal of radial length Sr, then 
the flame surface area fluctuation can be expressed as 


A ;(t) = E4 = x£4&(t)<5 r 


(23) 


where X k is the weight factors corresponding to the trapezoidal 
integration formula as given below 


A k - 


\ k = 1 and I< 
1 k # 1 or K 


(24) 


Linearizing Eq. (13) and taking Fourier transform leads to a first- 
order ODE for 'c(r) as: 


icaS(r) = u f -^ 


uj-V 2 f 


Further simplification reveals that 

T -ja)(r-Rfl) 

{(r) = ^ l - e v fV^ 

J CO 

Thus the fluctuating flame surface area is described as 


(25) 


(26) 


Af(m) = 


2Uf(co)n{R b -R a ) 


Ra 


S 2 t (co) 


(1 -cosS t (co) +j(sinS t (f«) -S t (a>))) 


S 2 ((o) 


(j (S t (co) cosS t (ro) - sinS t ((»)) - (1 - cosS t (co) 


-S t (ro)sinS t (oj))) 

where the flame Strouhal number S, is given as 

S t (co) = C ° {Rb ~ Ra) 


Vf\J~ Vf/uj) 


(27) 


(28) 


The flame transfer function T r (ft>) [14,27] between the unsteady 
heat release and the oncoming flow velocity can be shown as 


Tr(CO) = 


Q(oj)/Q A f (co)/A f 
Uf(co)/Uf Uf(a))/u f 


(29) 


where Af(co) and Af are given in Eqs. (27) and (21). Thus T r (S t ) can 
be rewritten as 


Tr(S t ) = 


2j 

(Rb + Ra 


S 2 


(S t - sinS t - j(cosS t - 1)) 


+ ^-(sinS f - S f cos S t - j(l - cos S t - S t sinS t )) 
S, 


(30) 


The characteristics of the flame transfer function T r (S t ) is discussed 
later in Section 4. 

2.3. Nonlinear self-sustained thermoacoustic oscillations 

Across the acoustically compact flame, the acoustic waves 
either side of the flame and the heat release fluctuation are related 
by mass, momentum and energy Eqs. (31)-(33). Substituting Eq. 

(31) into Eq. (33) and linearizing the momentum and energy Eqs. 

(32) and (33) give the time evolution of the outgoing traveling 
waves p u (t ) and pj(t) generated by the unsteady heat release Q'(t). 


Vr C Kb d H'lr t) C Kb a r z f 

A'(t) = . 1 / 2/ir^idr = X / £'(r, t)dr (22) ^nR 2 ’ pdz 

\ U f J R a U ' J Ra OL J z - 


= nRlipuU = 0 


(31) 
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8 ( 

dt 

8_ 

dt 


tzR 


nRl 


; f 1 pudz 

Jz r 

IK 


= nR h ( \p\Ij + pu[uti f ) = 0 


(32) 


pw yp 

2 y -1 


dz 


= nR 2 b 


y Z -i [P“]'/ +lpu[u 2 }^-Q(t)^ =0 
(33) 


By using these flow conservation equations across the flame and 
applying the end boundary conditions, the acoustic waves present 
in the combustion system can be solved. The resulting matrix equa¬ 
tion is given as 


A 


Pu(f)\ 

Pd + (f)y 




2(z f -zu) \ \ 
cud-Ml)J 
m-Zf) \ 


( 0 

Q'(t) 

\c u nR 2 b 


V 


(34) 


where A and B are coefficient matrices as given in Appendix B. Iter¬ 
atively solving the system equations enable time evolution of the 
flow disturbances to be calculated, thus providing a platform to gain 
insights on the onset of thermoacoustic oscillations as discussed 
later. 

Rayleigh index Z is generally used as a critical indicator of the 
combustion system stability. It describes the energy exchange rate 
between unsteady heat release and the pressure oscillations. If 
X > 0, the acoustical energy in the present system is increased 
and tends to drive the unstable process into a limit cycle. However, 
a negative Z indicates a ‘destruction’ interaction between unsteady 
heat release and pressure oscillations and so the combustion oscil¬ 
lations are decaying. Thus it would be interesting to check the Ray¬ 
leigh index Z. Following the work [15], we define a normalized 
Rayleigh index X(t) as given as 


Z (t) = J* d d(z-z f )dzj‘ 

p'(co, i//)Q'(co, iA)diA 


ill 


ill 


■2n/co «y2 


= di/^ 


Q. (ct>, i)dil/ 


(35) 


It is then used to characterize the heat-to-sound interaction in the 
present thermoacoustic system. 


3. Nonorthogonality analysis of thermoacoustic inodes 


Manipulating Eq. (36) leads to 


Pn 


-yp PVldp^ 

col pz \P 9z ) 


yp 


( 37 ) 


By using Eq. (37), the inner product between the eigenmodes p m 
and p* can be shown as 


r z ‘ 

(Pm,Pi) = / 

JZu 


PmPldZ = 


-yp f 

® 2 „ J*, 


z \ 9 nap; 


dz \p dz 


dz 


+ 


j( y-i) f Zi 

J Zu 


p m (T 8{z - z f )dz 


(38) 


To determine the inner product, unsteady heat release Q(co) needs 
to be estimated. If we assume the flame Strouhal number S, is small 
in the present system [17], then the classical time-lag J\f - z formu¬ 
lation can be used to describe the interaction between Q'(t) and 
oncoming flow [30], It has been demonstrated in Section 2 that 
the kinematic flame model can be linearized into the M - z formu¬ 
lation. In addition, this formulation has been widely used in quan¬ 
titative prediction of combustion instabilities [30], According to 
the M - z formulation, the unsteady heat release is given as 


Q'(t) 


ypUu' s (t - t v ) 


(y-i) 


(39) 


where A f denotes the interaction index describing the intensity of 
the interaction, u f (t) is the oncoming flow velocity at the flame 
position and r„ is the time lag between the velocity perturbation 
and heat release oscillation, which is assumed to be equal to 
0 Az d /u u . Here z d is the length of duct downstream of the flame 
holder. If we approximate u|(t- z„) in Eq. (39) to the first order 
i.e. u' f (t) - z v du' f (t) / dt, taking Fourier transform to both sides of 
Eq. (39) and substituting Q(co) into Eq. (38), then the inner product 
can be shown as 


(Pm,Pn) =Sr j' 

UJ n Jz u 


-yp r 

® 2 „ Jzr 


jAfc u e) <0nT ‘'Pj 
~ co n (R 2 b - R 2 a ) 


-jcom ( z-zj ) j aim ( z ~ z f ) 

P c u (l+Mu) -I- fp PQiO-Mu) 


?-C- d Mdz 

dz \p dz J 


h+ e l dO+ M d> 4- ni e l dO- M d) | — ()L (JEn 

Pd,m +Pd, m I dz\p dz 


Pus, 


Km) (pin 


-Pin 


dz 


(40) 


In order to gain insight on the stability behaviors of the present 
combustion system, nonorthogonality analysis of the combustion- 
excited modes is performed [25,26], The nonorthogonality is char¬ 
acterized by the inner product of two eigenmodes [25]. If the inner 
product is zero, then the eigenmodes are orthogonal. The combus¬ 
tion system with decaying orthogonal eigenmodes will be stable. 
However, when the inner product is non-zero, the combustion sys¬ 
tem is nonnormal and characterized with nonorthogonal eigen¬ 
modes. In such nonnormal system, finite-time growth of acoustic 
disturbances might trigger thermoacoustic instability. Thus it is 
critical to characterize the combustion system with or without 
nonorthogonal eigenmodes and to quantify the nonorthogonality. 

3.1. Nonorthogonal eigenmodes 


Let’s first consider the combustion system as described in 
Table A.l. The flame results in the mean temperature undergoing 
a jump, i.e. 7 d > T u . It behaves like an efficient monopole source 
[22 and produces the traveling acoustic waves as described in 
Eq. (5). Taking Fourier transform to it and taking the complex con¬ 
jugate yields 


-col 9 n dp;\ 

yp Pn dz \p dz j 


y —l (l*d(Z-Z f )(jCQ„)* 


(36) 


where co m and co„ are the mth and nth modes frequencies. They can 
be estimated by solving Eq. (34) in frequency domain as given by 
the resulting block matrix equation as 



/w, 

0 

0 

0 

0 ^ 


/ v ,\ 


/°\ 

C T 

0 

W 2 0 

0 

0 


v 2 

= c r 

0 

Si 

\ 0 

0 

0 

0 

VJ N ) 


\vj 


w 


w 


where V m = [p um (w m ),pj m (m m )] T ; c T = [1.1.... ,1]. Subscript m 
denotes the mth mode with frequency of co m and W m is a 2 x 2 
matrix as given as 


W,„ 


(42) 


where the time delay z u = 2(z f - z u )/c u (\ - M u 2 ) and 
z d = 2 (z d - z f )/c d ( 1 - M d 2 ) and the matrices A and B are given in 
Eq. (34). Theoretically, the acoustic disturbances consist of infinite 
modes. For simplicity, we assume that N modes are involved in 
the flow fluctuations, i.e. 
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N N 

Pu =XXm( a ’m), Pi = (43) 

m=l m=l 

Eq. (41) has solutions only if the determinant of the matrix W m 
vanishes, i.e. 

N 

det{W} = J] det{W m } = 0 (44) 

m=1 

The root satisfying det{W m } = 0, i.e. 

O m = co m +ja m (45) 

is the eigenfrequencies of the present system with the premixed 
flame confined. Physically, a m denotes the grow rate of thermoa¬ 
coustic oscillation. When a m > 0, the combustion system is unsta¬ 
ble. However, if a m < 0, the long-term evolution of eigenmodes 
are decaying. However, if the eigenmodes are nonorthogonal, they 
may interact and transient energy growth might occur. co m is the 
thermoacoustic oscillation mode frequency. They can be found by 
Newton’s iteration method [22], as shown later. 


4. Results and discussion 

The flame response to oncoming flow disturbances are studied 
first by using the model developed in Section 2. For a given forcing 
Uf(t), Eq. (13) can be numerically integrated in time to determine 
the evolution of the flame front. To illustrate the dynamic response 
of the premixed V-shaped flame, an unsteady inlet flow Uf(t) as 
given as Eq. (51) is applied to the flame. Here H(t) is the Heaviside 
step function and <f> 0 is the initial phase. 

1 +(1 +e)H(t)sin(cot+4> 0 ) (51) 

u f 

Fig. 2 shows the time evolution of the flame front. Comparison is 
then made with the mean flame position, as denoted by the dot line, 
it is interesting to observe that the flame moves upstream of the 
flame-holder during part of the oscillations as 2n/4 sg cut < 67t/4, 
since the fluctuating part of the oncoming flow velocity is larger 
than it’s mean value. 

As shown above, the flame dynamics can also be characterized 
by using the transfer function T r (S,). The variation of T r (S t ) with the 


3.2. Nonorthogonality arising from Robin boundary condition 


The forgoing analysis shows how to quantify the eigenmodes 
nonorthogonality of the choked combustion system. However, it 
is interesting to find out the contribution to the nonorthogonality 
arising from the complex boundary condition. This can be done 
by considering the combustor in the absence of flame, i.e. assuming 
AT —> 0. The nonorthogonality is also quantified by estimating the 
eigenmodes inner product (p m ,p„). Here p m and p n satisfy the 
homogeneous wave Eq. (5) in frequency domain as given as 


d 2 Pm 

dz 2 



(46) 


cPpn 

dz 2 



(47) 


By multiplying Eq. (46) with co n p„/ct> m and Eq. (47) with at m p m /<x> n , 
subtracting and integrating over the combustor, the inner product 
(Pm.Pn) is obtained namely 


/' 

Jz u 


d 2 p> 


d 2 p„ . 

dZ 2 P " ~ ()Z 2 P " 


dz = 


«“<) [ Zi 


f d 

/ fimPndz 

JZu 


Eq. (48) can be simplified to 

dpm . dp„ lZd 


dz Pn -fcPn, 


*-d nZ 

'■ + / 

. z u 7z u 

f Zd 

/ p m p n dZ 

J Zu 


d (dpm dpn dp n dpm \ j 

dz~dz~~dz~dz 


(48) 


(49) 


It is obvious that the first term involving the boundary conditions 
plays a dominant role in contributing to the inner product 
(Pm.Pn), while the second term on left hand side of Eq. (49) is zero. 
Thus the inner product (p m ,p n ) can be shown as 


r z d 

(Pm.Pn) = / PmPndz - 

JZu 


(®n-« 


dp. 


dpn , 


1 Zd 


dz Pn ~fc Pm 


(50) 


Eq. (50) can then be used to characterize the Robin-type boundary 
effect on the system nonorthogonality. The contributions to the 
nonorthogonality from the flame and the boundary are compared 
and discussed in the following section. 


(a) cot=0 (b) cot=7t/4 










Fig. 2. Time evolution of the premixed V-shaped flame surface anchored to the tip 
of the gutter during one period of the oncoming acoustic disturbance, as 
e = 0.2, R a /Rb = 0.25, co/2n = 59 Hz, V/ = 3.6 m/s, U/ = 37 m/s, and </> 0 = 1.7. 
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Fig. 3. Variation of flame transfer function T r with Strouhal number S,. 

flame Strouhal number S t is shown in Fig. 3. It can be seen that 
when the flame Strouhal number is negligible (i.e. S t —> 0 and 
m —> 0), the magnitude of T r is approximately 1 and the unsteady 
heat release and oncoming acoustic disturbances velocity is in 
phase, i.e. arctan{3(T r ), iR(T r )} = 0. However, as S, is increased, 
the amplitude of T r is decreased. And the phase between unsteady 
heat release and acoustic velocity is increased. Further increasing 
the flame Strouhal number S, gives rise to the amplitude of the 
flame transfer function to reduce to zero. However, the phase is 



Time (s) 


Fig. 4. Time evolution of (a) the phase diagram between Q(t)/Q and p'(t)/p; (b) 
Rayleigh index £(t). 


oscillating between -180° and 180°. The decrease of the magni¬ 
tude of T r with increased S t indicates that the flame is more sensi¬ 
tive to respond to lower frequency acoustic disturbances. However, 
high frequency acoustic disturbances pass through the flame more 
smoothly. 

The nonlinear response of the flame to the oncoming flow dis¬ 
turbances is then determined by numerically solving Eqs. (13), 
(16) and (34) with specified conditions as given in Table A.l. 
Fig. 4(a) shows the time evolution of the phase diagram between 
unsteady heat release Q(t)/Q and pressure fluctuation p'(t)/p, as 
denoted by the solid line. It can be seen that limit cycle is gener¬ 
ated at t = 0.5 s. The dash line describes the phase diagram 
between Q(t)/Q and p'(t)/p, which is assumed to consist of only 
the fundamental mode at co i. When the phase difference between 
them is n/2, circular-shape is observed. The mode growth rate is 
approximately e 6015t . 

Fig. 4(b) shows the variation of the Rayleigh index E(t) with 
time. It can be seen that when the initial perturbation grows into 
a limit cycle, Rayleigh index E is positive, keeps increased and 
finally ‘saturated’ at approximately 0.90. This verifies our numeri¬ 
cal model and explains why self-sustained combustion oscillations 
(i.e. thermoacoustic instability) occur. 

In order to quantify the eigenmodes nonorthogonality, the 
eigenfrequencies <o m need to be determined by using Eq. (44). 
Fig. 5 shows the variation of estimated coj with the flame location 
Z//L. Compared with the experimental measurements [29], qualita¬ 
tively good agreement is observed. It is worth noting that the mea¬ 
sured frequency at fixed flame location is changed with varied 
fuel-air ratio. 

With co m estimated, the inner product (p m ,p*) characterizing 
eigenmodes nonorthogonality can be calculated, by assuming that 
p- = 1.0 + jO, p+ = Z u p- exp{-jo)T u }, p+ = -p u W„/W 12 , and 
Pj = Z d p^ expj-jfuid}. Fig. 6(a) illustrates the variation of the 
inner product (p m ,p*) with the flame position z f /L. It is apparent 
that the thermoacoustic modes are nonorthogonal, since 
|| <Pm: Pn> II > 0. Furthermore, the nonorthogonality depends on 
the flame location Zf/L. It is also observed that the nonorthogonal¬ 
ity is increased with increased flame-acoustic interaction index AT, 
as shown in Fig. 6(b). This indicates that the system nonnormality 
becomes more pronounced with the interaction intensity between 
the flame and oncoming acoustic disturbances. 



Fig. 5. Comparison between the predicted mode frequency co, and the experimen¬ 
tally measured ones, as L = 1.92 m, M„ = 0.08, t„ = 0.01 s: - M = l.Oe - 5, 

— M — 4.5e - 4. 
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Fig. 6. The variation of || <p m , > ||/|| (Pi, Pi > II with the flame location z s /L, as 

t„ = 0.01 s, L = 1.92 m, and M„ = 0.08. c u = 340 m/m and T„ = 288 K. 


For the combustor as modeled in Section 2, the inlet is choked 
and so the mass flow rate is constant at z = z u , i.e. pu' + up' = 0. 
Further analysis shows that the acoustic impedance at this choked 
inlet Z(z u ) —p/u = -yp/u u ^ 0 and the boundary condition is a 
Robin-type one, i.e. Zdp/dz +j cop u p = 0. it is a weighted combina¬ 
tion of Dirichlet- and Neumann-type boundary condition [26], At 
the open outlet (corresponding to Dirichlet-type boundary condi¬ 
tion), the pressure fluctuation is zero, i.e. p(z d ) = 0 and so the 
acoustic impedance Z d = 0. The estimated inner product is shown 
in the bar chart as in Fig. 7. It can be seen that the acoustic modes 
are nonorthogonal. However, the order of magnitude of the inner 
product is much smaller than that in the presence of the flame. 
This indicates that the complex boundary condition of the ther¬ 
moacoustic system does contribute to the eigenmodes nonorthog¬ 
onality. However, the contribution is approximately 1.5% of that 
from the heat source—flame. In practice, there are many combus¬ 
tion systems such as baffled liquid-propellant rocket engines [1] 
or swirl-stabilized gas turbines, which involve with such choked 
boundary condition. These combustion systems are associated 
with nonorthogonal eigenmodes, whether or not a flame is present. 
This is different from those combustors with Dirichlet- or 
Neumann-type boundary conditions as analyzed in Appendix C; 


they are associated with orthogonal eigenmodes in the absence 
of the heat source. 

In general, eigenmodes nonorthogonality analysis reveals that 
the choked combustion system with the V-shaped premixed flame 
confined is a nonnormal thermoacoustic system. The nonnormality 
characterized by nonorthogonal eigenmodes arises from both the 
complex boundary condition and the flame. However, the contri¬ 
bution to the nonorthogonality from the heat source is identified 
to play a dominant role. The nonorthogonality is also found to 
depend on the flame location. In such a nonnormal system, small 
acoustic disturbances may exhibit transient growth before they 
eventually decay. If this transient growth is large enough, it can 
trigger combustion instabilities [31-46], 

5. Conclusions 

In this work, nonorthogonality analysis of a combustion system 
with a choked inlet and open exit is conducted. For this, a thermoa¬ 
coustic model is developed to study the nonnormal interaction 
between a premixed V-shaped turbulent flame and acoustic distur¬ 
bances. The state variables are expanded in terms of eigenmodes 
determined as the solutions of the linearized momentum and 
energy equations. It is shown the choked combustion system is 
nonnormal and characterized by nonorthogonal eigenmodes. The 
nonorthogonality is identified to arise from both the complex 
boundary condition and the heat source—flame. However, the con¬ 
tribution from the Robin-type boundary is approximately 1.5% of 
that from the flame. This indicates that the flame plays an important 
role in contributing to the system nonorthogonality. This is differ¬ 
ent from the combustion systems with Dirichlet- or Neumann-type 
boundary conditions, of which the nonorthogonality results mainly 
from the heat source. Such boundaries do not contribute to the 
system nonorthogonality. The eigenmode nonorthogonality results 
in the transient (finite-time) growth of acoustic disturbances, 
which may have great potential to trigger thermoacoustic instabil¬ 
ity. Classical linear stability theory is not enough to predict such 
finite-time system stability behaviors. Thus understanding and 
estimating the system nonorthogonality are important. 

Acknowledgements 

This work is supported by the Ministry of Eduction, Singapore 
under Contract No. AcR-Tierl-RG91/13-M4011228.050. This 
financial support is gratefully acknowledged. 

Appendix A. Geometry of the combustor 
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Fig. 7. The variation of || (Pm ,P„) || normalized with \\p m p n \ with m and n, as 
L = 1.92 m and M„ = 0.08. 


See Table A.l. 


Table A.1 

Geometry and flow parameters of the combustion system. 


Geometry 

Value 

Combustor length, L (m) 

1.92 

Radius of the combustor, R b (m) 

0.035 

Radius of the flame-holder/gutter, R a (m) 

0.0175 

Upstream length of the combustor, \Zf -z u \ (m) 

1.18 

Downstream length of the combustor \z d — z/|, (m) 

0.74 

Combustor inlet Mach number, M u 

0.08 

Combustor inlet temperature, T u (K) 

288 

Combustor inlet and outlet pressure, p u = p d (Pa) 

1.01325 x 10 5 

Fuel, 

Ethylene 

Mean heat release rate Q , (J) 

3.16 x 10 6 

Fuel-air ratio normalized on stoichiometric fuel-air ratio, 

0.7 

Flame speed (assumed constant) (m/s) 

3.6 

Combustion efficiency 

0.8 

Speed of sound at combustor inlet (m/s) 

340.17 
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Appendix B. Coefficient matrices 

The coefficient matrices A and B are given as 


* PdCd 

{^ + M 2 U - f(l-M u )(g-l) M u M d | + |l±* 


7 






(B.l) 




3-1 


id 




(B.2) 


Appendix C. Nonorthogonality analysis of combustors with 
Dirichlet- and Neumann-type boundary conditions 


Combustors are associated with two general types of boundary 
conditions, which are widely used in numerical and theoretical 
analysis. They are namely Dirichlet- or Neumann-type. In order 
to analyze the eigenmodes orthogonality of the combustors with 
such boundary conditions, Eq. (50) can be used: 


(1) Dirichlet (corresponding to open end), i.e. p(z u ) = p(z d ) = 0. 
The inner product (p m ,p„) as indicator of orthogonality can 
be shown as 


r Zi 

(Pm,Pn) = / PmPndZ 
Jz u 



dp m _<9p„ 
dz Pn dz Pm 


2d 

Zu 


= 0 


(C.l) 


It can be seen that the combustors with the Dirichlet-type 
boundary conditions are associated with orthogonal acoustic 
modes. The acoustic impedance Z = p'/u' is zero at the open 
end in the absence of heat source. However, when a heat 
source is present, the combustion systems become nonnor¬ 
mal characterized with nonorthogonal eigenmodes [25]. This 
indicates that nonnormality arises from the heat source, 
instead of the boundary conditions. 

(2) Neumann-type boundary condition, i.e. u'(z u ,t) = 0 and 
u' (z d , t) = 0. 

Analyzing the nonorthogonality of combustors with 
Neumann-type boundary conditions begins with estimating 
the inner product of the eigenmodes (p m ,Pn) as given in Eq. 
(C.2). It can be shown that (p m ,p n ) is zero. This means that 
orthogonal modes are involved and the system is normal. 
The acoustic impedance at the close end is infinitely large 
Z — * oo. 
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